Coulomb correlation in presence of spin-orbit coupling: application to plutonium 
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Attempts to go beyond the local density approximation (LDA) of Density Functional Theory 
(DFT) have been increasingly based on the incorporation of more realistic Coulomb interactions. 
In their earliest implementations, methods like LDA+[7, LDA + DMFT (Dynamical Mean Field 
Theory), and LDA+Gutzwiller used a simple model interaction U. In this article we generalize the 
solution of the full Coulomb matrix involving F' ' to F' 6 ' parameters, which is usually presented 
in terms of an Imt basis, into a jmj basis of the total angular momentum, where we also include 
spin-orbit coupling; this type of theory is needed for a reliable description of /-state elements like 
plutonium, which we use as an example of our theory. Close attention will be paid to spin-flip 
terms, which are important in multiplet theory but that have been usually neglected in these kinds 
of studies. We find that, in a density-density approximation, the jmj basis results provide a very 
good approximation to the full Coulomb matrix result, in contrast to the much less accurate results 
for the more conventional £me basis. 

PACS numbers: 71.10.+X, 71.20.Cf, 64.60.Cn 



I. INTRODUCTION 

Strongly correlated electron systems are solids where 
the important outer-shell electrons have two conflicting 
and opposite tendencies. On one hand, they maintain a 
strong memory of the atomic or localized orbitals from 
which they arise, which have a large electron-electron 
electrostatic interaction between discrete states. On the 
other hand, the same electrons hybridize with neigh- 
boring orbitals causing them to delocalize by tunnelling 
from one atom to those nearby, forming chemical bonds 
and spreading out the discrete atomic states into nar- 
row energy bands. In such systems, it is therefore neces- 
sary for a correct description and understanding of their 
electronic properties to maintain both of these aspects. 
The second tendency is very accurately calculated by 
density functional theory (DFT) band-structure calcu- 
lations, while the first one involves a consideration of 
many-body effects, and, while more difficult to treat, is 
still a crucial aspect of the physics. Thus, it is important 
to increase our knowledge of the Coulomb interactions 
that strongly affect the atomic character of these sys- 
tems, an effect which is often underestimated or poorly 
approximated in calculations that include details of the 
band-structure. 

These types of effects are particularly important for 
the electronic structure of /-electron elements in gen- 
eral, and especially for the actinides. For these materi- 
als, density functional theory (DFT) calculations in the 
local density approximation (LDA)i^ often give signif- 
icant discrepancies. For example, 5-Pu from this kind 
of approach is predicted to have an equilibrium volume 
25% smaller than experiment, which is the largest known 
deviation from LDA. To overcome these difficulties, var- 



ious attempts to go beyond LDA have been proposed, 
such as I. DA • / LDA+DMFT,£&I and more recently 
LDA+Gutzwiller J£ All three of these methods add a lo- 
cal Hubbard-likc term to a band Hamiltonian, and re- 
quire subtraction of an average LDA Coulomb interac- 
tion (the double counting correction). The differences 
between the various methods reside in the way the ef- 
fects of this interaction term is handled. 

In the LDA+J7 approach, which employs a Hartree- 
Fock mean-field solution, the Hubbard term leads to an 
orbital-dependent shift in the potential. Such a crude 
mean-field approximation is questionable for cases involv- 
ing strong correlations. In the LDA+Gutzwiller method, 
a variational wavefunction is built, for which the mean 
values of the interaction are calculated exactly. In the 
more sophisticated DMFT method, the effect of the in- 
teraction is described by a self-energy, which acts as an 
energy-dependent complex potential. In this approxi- 
mation the self-energy is assumed to be local (i.e., mo- 
mentum independent) and is determined self-consistcntly 
within an impurity-like approach of a correlated site em- 
bedded in a effective bath. In many DMFT calculations 
the Hubbard-like term, at least in the early implementa- 
tion of these methods, has often been treated in a fairly 
simplistic way. For example, some applications use a sin- 
gle U term, average over all interactions, while others 
also include an exchange-averaged parameter J. Over 
time the general tendency for all three methods has been 
to include more and more realistic interactions. For ex- 
ample, since its first use, the LDA+C7 has usually been 
rotationally invariant^ However, a multiband version of 
Gutzwiller approach developped by Bunemann, Gebhard 
and Weber— was able to handle spin-flip terms and very 
recently, a multiband generalization of slave-boson for- 
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malisroii has been proposed to be rotationally invari- 
ant^ To make progress, it is clearly important to de- 
velop a more sophisticated treatment of the Coulomb 
interactions. In addition, for high-Z materials like the 
/-electron actinides, spin-orbit must also be accurately 
included. To do this is the goal of this article. 

The paper is organized as follows. Section II is devoted 
to the presentation of Coulomb matrix elements in the 
jrrij basis. In this section, we first formulate a general 
expression for the interactions. We then show how to 
make an approximate density-density correlation calcu- 
lation for these interactions. The corresponding matrix 
elements will be tabulated in terms of Slater integrals. 
Section III presents the eigen-spectrum of the atomic 
Hamiltonian in various different approximations. In Sec- 
tion IV we use these eigenvalues to study the single par- 
ticle spectral density. In particular, the quality of various 
approximations will be evaluated against a rigorous solu- 
tion. Finally, we conclude in Section V by stressing the 
main results of our approach. 



II. COULOMB MATRIX ELEMENTS IN THE 
jm 3 BASIS 



the Schrodingcr formulation as a limit of the Dirac equa- 
tion. It is due to the interaction of the magnetic moment 
of electron spin with the effective magnetic field created 
by the orbital motion, and has the following expression: 

N 

H so = J2^)LA , (2.3) 

with, by dropping the index z, 

1 dV 

« r > = 2^*7F- (2 - 4) 

The spin-orbit interaction is diagonal in a jrrij basis and 
splits the j = I ± \ states into two subsets having correc- 
tion energies A £+ i = £xnt/2 and A £ _i = — (£+l)x n e/2, 
respectively. The splitting energy is (21 + l)xnt/2- Here 
Xn£ is given by the radial integral: 

Xm = fr 2 j \{r)\R nl {r)\ 2 dr . (2.5) 

Spin-orbit coupling begins to be important for atoms 
with atomic number Z > 20 , where the derivative Up- 
starts to become significant. 



A. Electronic structure of an isolated atom 

The Hamiltonian H of an isolated many-electron atom 
or ion is 

r-f \ 2m n J ~ In- rj\ 

i= i x / %>j 1 j 1 

(2.1) 

Beside the interaction of the electrons with the posi- 
tive charge (Z) of the nucleus, this Hamiltonian contains 
two important features: the spin-orbit coupling and the 
electrostatic (Coulomb) interaction between electrons. 
By neglecting spin-orbit, in the central field approxima- 
tion,— the eigenstates of the system are Slater determi- 
nants built from individual states \n£mesm s ) having the 
following wave function: 

<t>nt mismB {r, 6, ip) = ^ r Yi mi (9, <p)r] sms (s z ) . (2.2) 
r 

Here Yi me is a spherical harmonics, R n g the solution of a 
radial Schrodingcr equation, and rj smg (s z ) an eigenfunc- 
tion of S z . The set of indices nlmgsm s is sufficient to 
determine completely a state with eigenvalue E n g, having 
the degeneracy 2(2£+l). It is this basis (or its equivalent 
in a jrrij basis) that we will use to study the full Hamil- 
tonian (|2.1[) . The appearance of a two-body term, i.e., 
the electrostatic interaction between electrons, makes the 
problem sufficiently complicated so that an cigenstatc, 
even in a pcrturbative description, will not be in general 
be a single Slater determinant. 

The spin-orbit coupling, which is still a one-body oper- 
ator, is a relativistic effect and can be directly obtained in 



B. Two limiting behaviors: LS or jj couplings 

Depending on which of two contributions, the Coulomb 
or spin-orbit interaction, dominates, there are two lim- 
iting regimes. The first is the LS coupling or Russell- 
Saunders regime, in which the electrostatic exchange 
interaction is predominant; this is responsible for the 
Hund's rule ordering of states. In this case Imt is the 
most convenient basis with the unperturbed states in the 
form \n£m£sm s ), since the Coulomb matrices are diago- 
nal in spin. 

The opposite regime, the jj coupling regime, occurs 
when spin-orbit coupling splitting is greater than the 
electrostatic terms. In that case, it is convenient to work 
in the jrrij basis, which diagonalizes the spin-orbit term 
and to treat the electrostatic term as a first-order per- 
tubation, leading to a single (diagonal) correction to the 
unperturbed cigenenergics. 

For actinides we are in an intermediate regime. For 
example, the average exchange J in plutonium^i is of the 
order of 0.7eV and the spin-orbit parameter X5f f° r 5/ 
states is in the range 0.25-0.54 eV, producing an energy 
splitting between j = 5/2 and j = 7/2 states in the range 
of 0.9-1.95 eV. In this case, it is important to treat the 
spin-orbit and the electrostatic terms on the same footing 
by diagonalizing them in a given basis. Transformation 
from one basis to the other can be performed with the use 
of the Clcbsch-Gordan coefficients (£ me s m s \j rrij) with 

li m i) = \^ m e s m s ) (£ m e s m s \j rrij} . (2.6) 

As we will explain below, there are strong arguments for 
using the the jrrij basis. In this case, the diagonal part 
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gives directly the jj coupling approximation and the fully 
diagonalized result provides a reliable description of the 
intermediate regime. 

C. Coulomb interaction in the jrrij basis 

In the jrrij basis we can write the Coulomb contribu- 
tion to the Hamiltonian as 

V C = \ Yl ^1234C t 1 ct C4 C3 (2.7) 
1234 

Here 1, 2, 3, and 4 are a shorthand for individual par- 
ticle states \n£ jrrij). The spatial part of the Coulomb 
interaction can be expanded in the £mi basis as 

1 °° fe r k Air 

1 1 k=0q=-k > 

Here r< (r>) is the lesser (greater) of r and r', and Yk,q 
is a spherical harmonics, with the solid angle spanned by 
ft = (6,<p). 

The Coulomb matrix clement V1234 is explicitly given 
by (using the system of units where e 2 = 1) 

^1234 = (jirn 1 j 2 m2\-, — — T , \ j^rn^jirrn) . (2.9) 
\r — r I 

In this expression different levels of approximation can be 
made. If it is used exactly as is with no approximation, 
we will refer to the results as involving "spin- flip" terms, 
since the creation and destruction operator in Vc can 
flip spins. We will call the next level of approximation 
the "density-density correlation" approximation since, as 
it will become clear below, in this approximation one 
retains only the case for which either 1 = 3 and 2 = 4 
(the direct term) or 1 = 4 and 2 = 3 (the exchange term) . 
Thus this part of the Hamiltonian reduces to 

Vc = \ ^(^1212 -V 1221 )mn 2 (2.10) 

1234 

It is worth noting that the usual selection rule that occurs 
in the £rri£ basis, namely, that the exchange interaction 
vanishes for antiparallcl spins, does not occur here, since 
the jrrij basis is a mixture of different luii and sm s states. 
As a result, in our present case, a net interaction within 
the density-density correlation approximation is always 
the difference between a direct and an exchange term. 

The density-density correlation approximation is very 
important because it make the LDA+U feasible. For 
the case of DMFT, if one uses the Hirsch-Fye typei& 
Quantum Monte Carlo (QMC) solver, the Hubbard- 
Stratonovitch transformation makes this approximation 
necessary. For the Gutzwiller case, even if it is in prin- 
ciple possible to keep the spin-flip terms^ it is however 
much easier to avoid these terms and to make the density- 
density correlation approximation, especially for the 
density-matrix derivation of the generalized Gutzwiller 
method^ 



1. General formulation including spin-flip terms 

Our starting point is the definition (|2.9[) . Since we 
are mainly interested in /-electron elements (£ = 3), the 
possible values of the j = £± | in this expression is either 
7/2 or 5/2. In the Coulomb potential expansion (|2.8p . it 
is suitable to insert a closure relation in the decoupled 
basis, and make further use of selection rules of Clebsch- 
Gordan coefficients. As a result, we obtain (see Appendix 
for a detailed demonstration): 

(jim 1 j2m2\V\j 3 m 3 j4m 4 ) = (-) m3 ~ m2 <5„ ll+m2:I „, + „ l4 

xt F(k) ^§ <^0|fc0) 2 

x Blijxmi, j 3 m 3 )B e k (j 2 m 2 ; jim^) , (2.11) 

where the summation over k extends only over even val- 
ues of k ranging from to 21 due to selection rule for the 
Clebsch-Gordan coefficient (£0£0\k0) (see Appendix). 
The Bf.(jm; j'm') are given by 

Bi(jm;j'm') = £ 2a(£ (m-a) l - a\jm) 
x (£ (m'—a) — a\j'm') 

x (£ -(m-a) £ (m'-a)\k {m'-m)) , (2.12) 
and the F k are the Slater integrals 

= f drdr'\R nl {r)\ 2 ^\R n( {r')\ 2 . (2.13) 
J r> 



2. Density-density correlation approximation 

As explained above, in the density-density correlation 
we retain among all possible terms of (|2. 1 1|) only those 
for which (j'iTOi) = (j 3 m 3 ) and (j 2 m 2 ) = (.747714) for 
the direct term, whereas (j\mi) = (J47714) and (j 2 m 2 ) = 
(J37713) gives the exchange term. By using (|2.1ip with 
these restrictions, we obtain these two kinds of matrix 
elements: 

^mfms = {jimi j 2 m 2 \V\jimi j 2 m 2 ) , 

J mlmi = {jimi j 2 m 2 \V\j 2 m 2 jimi) , (2.14) 

which are found to be 

21 

k=0 
21 

J tt 2 m2 = E^OimijjamaJirf*), (2.15) 

k=0 
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where 



(Of _i_ i ^2 

4(hm i; hm 2 ) = (-)^-^ L^_ (^0^0|fcO) 2 



x B i (h m i i h m i ) B l ( h m 2 ; h m 2 ) , 

(2.16) 



and 



(21 4- l) 2 

6i(j 1 m 1 ;j 2 m 2 )= (2fc+1)2 (^0^|fc0) 2 

xBi(j imi ; j 2 m 2 ) 2 . (2.17) 

We have checked that we can use these two formulas 
for £ < 2 to retrieve results first established by Inglisi^ in 
the early 1930's and which are now common in textbooks 
(see, for example, Ref. [Til ). Since results for /-electron 
elements (I = 3) in these references are not provided, we 
have tabulated them in Tables [IT] and [HT] of the Appendix, 
where we give also some symmetry relations they obey. 
These results are a generalization for the jrrij basis of the 
more familiar imc results: 



21 

fe=0 
2i 

k=0 



(2.18) 



with 



4tt 



bk 



+ T E (Y em \Y k * q \Y em ){Y em ,\Y kq \Y em ,) 
Air 

7 

q— — k 



q— — k 
+ k 



2k 



■I 



(2.19) 



3. Averaged interactions and the single U limit 

The direct calculation of F 1 -^ from the atomic wave 
functions in their definition (|2.13j) usually overestimates 
their values, since it neglects screening effects that occur 
in the solid, i.e., the relaxation of other electrons pushed 
away by the two interacting electrons, which leaves be- 
hind a net positive charge and reduces the strength 
of the interaction. Consequently, effective values for 
F 1 -^ arc usually computed from constrained LSDA (lo- 
cal spin-density approximation) DFT calculations.— Un- 
fortunately, because LSDA can not distinguish between 
individual orbitals, it is only sensitive to the total spin- 
density and hence this approach can only provide an av- 
eraged direct and exchange interactions. 

For example, if this average is performed for the £mi 
basis, one obtains for /-states the following averaged di- 



rect and exchange interactions^ 



u = 



(2£+l) 2 
F (o) 



(2.20) 



and 
U-J 



>ZTT)^>" 



J mm' ) 



2£(2^+ 1) 

= F (0) - (286F (2) + 195F (4) + 250F (6) )/6435 . 

(2.21) 

^From the constrained LSDA values of U (= F( ') and J, 
and requiring constant ratios of /F^ and F^ /F^ , 
it is possible to assign a unique value to each F^ (k=0, 
2, 4, 6). This procedure, when applied to plutonium, 
leads to the values given in Table [1^ 

We have generalized the same kind of average for the 
direct and exchange interaction in the jrrij basis for / 
states, and find: 



U = 



1 



(2j 1 + l)(2j 2 + l) 



£« =i?(0) < (2-22) 



and 



J 



3132 



(2j 1 + l-S hj2 )(2j 2 + l) 



7j£«- (2-23) 



The exchange terms vary if they are taken for different 
h~h Pahs (5/2-5/2, 5/2-7/2, and 7/2-7/2). Accord- 
ingly, one finds: 



^5/2 5/2 = 


30^35 


+ ^ (4) ), 




•^5/2 7/2 = 


J7/2 5/2 = 


! ( A F (2) 

48 v 25 


40 
+ 77- 


J 7 / 2 7/2 = 


1(12 F (2) 

56 v 21 


+ h 
77 


200 
" 429 



143 



(2.24) 



If we go one step further in the averaging process and 
use a single J value, regardless of which j (5/2 or 7/2) it 
comes from, then we obtain: 



J = J_(^ F (2) 

182^32 



28^ + ^(0). (2.25) 
11 429 ' 1 ' 

These values can be used to further approximate the 
density-density approach, since the obvious simplification 
is to replace the detailed interaction given by (|2 . 14[) by a 
single value, which is their average over all possible pairs 
of orbitals. This last approximation has been widely used 
and we discuss next its effect on the spectrum. 



III. SPECTRUM OF EIGENVALUES 

A. Hamiltonian and Fock states 

We now explain how to obtain the eigenvalues of the 
atomic local part of the Hamiltonian (i.e., with no kinetic 



TABLE I: Slater integral values (in eV) for Plutonium from 
Ref . H 



^(0) 


i?(2) 


pW 


J?(6) 


4. 


8.343639 


5.57482 


4.12446 



energy term) in second quantization form: 

Hioc = - ^2 ^I234ct4c 4 c 3 Aicjci , (3.1) 



12.1 1 



where p (=1, 2, 3, or 4) is a shorthand for j p m p , V1234 is 
the Coulomb potential (see Eq. 12. 9p , and Ai is the spin- 
orbit term (see just below Eq. 12. 4j) . which is diagonal in 
the jrrij basis. 

First of all, this Hamiltonian conserves the number of 
particles. Thus, for a given occupancy N (0 < N < 14) of 
/-levels, the dimension of the basis for configuration f N 
is 14!/iV!(14 — iV)!, which is the number of Fock states 
njy) for spreading N electrons among the i = 1 to 14 
individual quantum states arising from the six 5/2 states 
(with rrij ranging from -5/2 to +5/2) and the eight 7/2 
states (with rrij ranging from -7/2 to +7/2): 



M = (c t 1 ) 0l ...(4) 014 |o) 



(3.2) 



Each state of the form (|3 - 2[) in this subspace has an 
electron occupancy of N, i.e., J2°i = N, where Oi is ei- 
ther or 1 in order to obey the Pauli principle. We com- 
pute then the matrix elements (n' N \H\n^) of the Hamil- 
tonian in this Fock basis of states. The diagonalization 
of this matrix gives the energies Ea of the multiplet 
structure and their corresponding eigenstates \a), with 
the Fock state components: 



\a) = ^(nAr|a)|nAr) . 



(3.3) 



As mentioned above, if only the diagonal part of Hamil- 
tonian, which corresponds to density-density approxima- 
tion, is kept (valid in the limit of the jj coupling), the 
Fock states are already eigenstates of the system, and the 
eigenvalues are the diagonal matrix elements. 

Moreover, as we are dealing with isolated atoms, the 
eigenvalues are discrete. Thus, to plot the spectrum, we 
have added to each eigenvalue a small imaginary part, 
which modifies the spectra from a sum of d functions 
to a sum of equal-width Lorentzians (here set equal to 
0.1 eV). 



B. Results 

In Figs. [1] and the density of states (DOS) is the sum 
of all Lorentzians from each eigenvalue. We note that this 
DOS should not be confused with the spectral density 
presented in next part: it is just a visual way of present- 
ing the distribution of eigenvalues. For illustration, we 
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FIG. 1: (Color online) Atomic spectrum of eigenvalues of 
Pu without spin-orbit coupling, with N=5 electrons for single 
U (dotted line), density-density approximation (dashed line) 
and spin-flip terms (full line) included respectively. 



present the case of plutonium where the Slater integrals 
parameters are given in Table I. Figure [1] presents the 
atomic spectrum of eigenvalues of Pu without spin-orbit 
coupling, for a given occupancy of N = 5 electrons. The 
different approximations with increasing order of com- 
plexity are displayed: single U , density-density approxi- 
mation and spin-flip terms included respectively. Figure 
[2] displays the same situation when spin-orbit coupling 
has been taken into account for occupancies of N = 5 
and N = 6 electrons. 

The most accurate description is given when the spin- 
flip terms are not neglected. In this case, the spec- 
trum is much more structured than the result for the 
density-density approximation. The spin-flip terms in 
the Hamiltonian couple different states that arc degener- 
ate in the density-density approximation, and split them. 
This causes the additional structure in the spectrum. 
The density-density approximation spectrum therefore 
has less structure, but has its gross features centered 
in the same region as the spectrum which include the 
spin-flips. For the single U limit result, for a given occu- 
pancy N and neglecting of spin-orbit coupling, there is a 
single eigenvalue UN{N — l)/2 whose degeneracy equals 
the dimension of subspace. This is the origin, when 
smeared with a Lorentzian, of the single peak, which is 
actually too high in energy and well separated from the 
more realistic spectrum in the density-density approx- 
imation or the one including spin-flip terms (roughly 5 
eV above the center of gravity of the other ones). Clearly, 
the single U limit overestimates the interaction energy. 
This is mainly caused by the complete neglect of any ex- 
change term, which would reduce the average interaction. 
This approximation, which has been widely used, can be 
a problem for an atom embedded in a solid, since the 
metal-insulator transition is a delicate balance between 
localization due to the interaction U and derealization, 
represented by a bandwidth W. Thus the single U limit 
could artificially push the system towards an insulating 
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FIG. 2: (Color online) Atomic spectrum of eigenvalues of 
Pu with spin-orbit coupling, with N=5 and N=6 electrons 
for density-density approximation (dashed line) and spin-flip 
terms included (full line) respectively. 



state, or at least, increases the localized characters of the 
/-electrons. To remedy this situation, without increas- 
ing the level of complexity of the solution, would be to 
replace the single U limit by a single U — J, since we have 
seen above that there is always an exchange term in the 
density-density approximation for the jrrij basis. In that 
case, taking the overall average exchange J = 0.323eV 
of Eq. (|2.25p . one would obtain a single peak located at 
36.77 eV, which is in better agreement with more elab- 
orated results, and especially very close to the central 
maximum of the density-density result. 



IV. ATOMIC TEMPERATURE GREEN'S 
FUNCTION 

As mentioned in the introduction, one cannot neglect 
a correct description of the atomic aspects of strongly 
correlated electrons systems. In this part we concentrate 
on the atomic (i.e., local) Green's function of an atom, 
which will be later embedded in a solid. The results 
presented in this section can be seen as a first step, with 
further hybridization to the rest of the medium to be 
be added later as, for example, is done in the DMFT 
approach. We start from the definition of (imaginary 
time) temperature Green function : 



9*(t) = -<T rCfc (T)4(0)> , 



(4-1) 



where ct is a creation operator in a state k, which is, in 
the present case, one of the 14 atomic orbitals of the jrrij 
basis for I = 3. Based on the antiperiodic property 



S fe (T + /3) = -S fc (r) , 



(4.2) 



where 1 / (3 is the temperature, one can write the Fourier 
expansion of Bh(j) as 



Sfe(r) 



r S(iw n ) , 



(4.3) 



with the Fourier transform 

%{iu n ) = / e* 
Jo 



r Sk(T)dT 



(4.4) 



The Matsubara frequencies u n are given by: 

u n = (2n + 1)tt/P . (4.5) 

If the cigenstates \a) and their eigenvalue E a are known 
from a diagonalization of the Hamiltonian, then Eq. (|4.4|) 
has the following Lehmann representation: 



1 



-PK 



a\Ck\ 



EL 



Here Zq is the grand partition function: 



(4.6) 



(4.7) 



E' a = E a — (J,N a , and /i, the chemical potential, is chosen 
to fix the average number of particles with: 



(4.8) 



Eq. (|4.6[) can be used to define the spectral density; this 
can be compared to photoemission experiments, and is 
given by 



J-^|(a| Cfe |6)| 2 (e-^+e-^)5(a;+K-^) 



a . b 



(4.9) 

Because the sum involves discrete cigenstates for the 
atomic case, we add to oj a small imaginary part so that 
the spectral density becomes a sum of Lorentzians when 
the results arc presented. It should also be noted that 
for solids the atomic Green's function tends to become a 
very good approximation at sufficiently high temperature 
(see, for example, the discussion of Fig. 2 in Ref . [l9t) . 



A. Comparison of the single U model with the 
density-density approximation 

Since all states of a given occupancy N are degenerate 
with eigenenergy UN(N — l)/2 when spin-orbit is ne- 
glected, the single U model has an analytical expression 
for the Green function, which can be written: 



-E 



13! 



Z G ^ JV!(13 - N)\ 

-PN^U-A , -0(N+l)[%U-n] 
X " : ,rrr ■ ( 4 '10) 

This result, which can be presented in a form previously 
obtained by McMahan et al.^ shows evenly spaced val- 
ues separated by U in the poles of the Green function. 



10 



O 4 
Q 




-10 



-5 





ENERGY (eV) 



10 



FIG. 3: (Color online) Spectral density of atomic Pu at 
T=15800K, for single U limit (full line) and density-density 
approximation (dashed line) respectively, including spin-orbit 
coupling in both cases. 



When spin-orbit is switched on / 0), it is less con- 
venient to give an simple analytical formula, since the 
degeneracy of the states is partially lifted. We present 
the spectral density for this case in Fig. [3] One can 
still see the effect of approximately evenly spaced values 
of the order of U for the poles, but additional structure 
also appears. Note that the effect of temperature, here 
chosen to be 15800 K, allows additional excited states 
to appear in the spectra. A lower temperature would 
have completely quenched much of this additional struc- 
ture. In addition, this temperature has been used in some 
recent DMFT Pu calculations to check the Quantum 
Monte Carlo (QMC) DMFT codes against atomic-like 
Hubbard-I approximations^, since the atomic limit is a 
good approximation at this temperature. The average /- 
occupancy is chosen to be about 5.6 in accordance with 
the result of Re&. The comparison with the density- 
density approximation spectral density (dashed line) in 
Fig. [3] clearly shows that the single U approximation is 
an oversimplified model for the correlation effects, and 
misses many of the features of even the density-density 
results, and misplaces some of the peaks. 
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FIG. 4: (Color online) Exact results, which include spin-flip 
terms in Hamiltonian, for the total (full line), and 5/2 (dot- 
ted line) and 7/2 (dashed line) projected spectral densities of 
atomic Pu at T=15800K, including spin-flip terms in Hamil- 
tonian 

In Fig. [5] we compare the spectral density for the 
density-density approximation with the exact results for 
both j-j and Is coupling. Both show a similar overall 
gross structure of the spectrum, which is located in ap- 
proximately the same energy range. However, the j-j 
coupling clearly does a much superior job in reproducing 
the actual peak structure in the spectrum. In addition, 
the excellent agreement between the j-j density-density 
approximation and exact spectral densities suggests that 
in many situations the j-j density-density approxima- 
tion may be an adequate approximation to more refined 
treatments of the Coulomb effects for Pu. It is worth 
mentioning that, when all terms are retained as here, i.e., 
spin-orbit and full electrostatic interaction, the results do 
not depend on the basis chosen for the representation; the 
local Hamiltonian that would have been solved from the 
LS basis would have had the same final spectral density. 
Thus the similarity we stress between the full result and 
the density-density approximation is not an artifact of 
the choice of representation, namely the jrrij basis. 



B. Comparison of the density-density 
approximation with exact results 

In Fig. 2] we show the exact results for the spectral 
density (within the basis we have chosen). These include 
the spin-flip terms and involve exact diagonalization of 
the local atomic Hamiltonian. In this figure we also in- 
clude the decomposition of the spectra into its 5/2 and 
7/2 contributions, with their respective weights of 6 and 
8. As expected, the spin-orbit splitting has displaced the 
gross features of the 5/2 density downward with respect 
to those of the 7/2 features. These spectral densities have 
much more structure in them than the density-density re- 



V. SUMMARY AND CONCLUSION 

We have presented a study of the atomic electronic 
structure of the /-electron elements in the presence of 
spin-orbit coupling, and applied this to plutonium. We 
have emphasized the role of electronic interactions in a 
jrrij basis. For this purpose, we have derived a general 
analytical expression for the matrix elements of the inter- 
action in this basis; some of these that can be used in a 
density-density approximation have been tabulated, since 
they could be useful for possible calculations in solids. 
We have used these to study the effect of Coulomb in- 
teractions on the eigenvalue and spectral functions. Di- 
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FIG. 5: (Color online) Comparison of the spectral density of 
atomic Pu, including spin- flip terms in Hamiltonian (full line), 
with the density-density approximation result (dashed line). 
The upper panel (a) provides results for the density-density 
approximation for j-j coupling (rijUji) whereas the bottom 
panel (b) is for the case of Is coupling (n; s n;/ s /). Note that 
for the latter case the full matrix must be diagonalized since 
the spin-orbit coupling term is not diagonal in the Is basis. 
The spin-flip case is reproduced in both panels for comparison. 



agonalizations of this Hamiltonian for different occupan- 
cies and computations of spectral density, have been per- 
formed in various approximations, including single U, 
single U — J, density-density ,and, finally, the full interac- 
tion matrix retaining all spin- flip terms. We have found 
that the density-density approximation, in the jrrij basis, 
which neglects off-diagonal terms (i.e., spin-flip terms) in 
the interaction matrix, gives excellent agreement with the 
full interaction. 

The high accuracy of the j-j density-density approxi- 
mation should make it very useful in electronic-structure 
calculations for solids. For example, it is a very tractable 
method for the Hirsch-Fyc algorithm of the QMC solver 
in DMFT, and should makes the application of the 
Gutzwiller method much easier. 
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of three spherical harmonics: 



APPENDIX: EXPRESSION OF THE GENERAL 
MATRIX ELEMENT IN THE jrrij BASIS 



We demonstrate now formula (|2.1ip . The general ex- 
pression element of the Coulomb interaction in jrrij basis 
is V1234 = (jiMijiM^VljsMsjtiMi). In this appendix to 
avoid too heavy notations, we will use the notation \ jM) 
for an element of the jrrij basis and \£ma) to designate 
the state \£mi sm s ) of the £mi basis. j\ to j'4 could be 
one of either £ + \ or I — | total momentum quantum 
number (£ being well defined here). By inserting four 
times the Clebsch-Gordan (CG) expansion: 



\jM) = \£ma)(£ma\jM) 



(A.l) 



in the definition of V1234, where (£ma\jM) is a shorthand 
for the (chosen real) CG coefficient (£m^a\(£^)jM), one 
obtains: 

V1234 = 2J (£mi<Ti\jiM 1 )(£m2<T2\j2M 2 ) 

x (£m 3 a 3 IJ3-M3) (£m 4 a 4 \j4M4) 
x(£miai £m,202\V\£m 3 a 3 £7714(74) . 

(A.2) 

The matrix element of interaction in the £mi basis ap- 
pearing in this last formula can be expressed as: 

{£mioi£m2(J2\V\£rn3<j3£rn i <7 i ) = 5a 1 a 3 5 l j 2 a l 

x J drdndr'dn'Y; mi (n)Y; m2 (n') 

x \R nl (r)\ 2 —^- 



Tl \Rndr')\ 2 Y em3 (n)Y em4 (n') 



(A.3) 

^From Eq. (|2.8[) . the above integral can be conveniently 
cast into a sum of a product of one radial integral (iden- 
tified as the Slater integral (|2.13[) ) and two angular inte- 
grals: 



00 fe 4-7T C r k 
EE^TI drdr'\R ne (r)\^\R ne (r')\ 



x J dnY e * mi (n)Y k * q (Q)Y im3 (Q) 

x f dn'Y e * m2 (n')Y kq (n')Y emi (n') . (A.4) 



With conjugation relation Y^ m (fl) = (— ) m Y e _ m (fl), 
each one of the two angular integrals can calculated from 
a well-known identity giving the integral of the product 



dQYe imi (il)Yi 2m2 (il)Yg 3m3 (£l) — (— ) 



)>!;; 



'(2£x + 1)(2£ 2 + 1) 



V 4ir(2£ 3 + 1) 
x (40^0140) (£ 1 m 1 £ 2 m2\£ 3 -ma) 



(A.5) 



Application of this last identity with the selection rule 
for CG coefficients (^ittii^^Ks"^) which vanish unless 
7773 = mi + rn 2 , enables to rewrite (|A.3[) : 



(£m 1 a 1 £m2a2\V\£m 3 a3 £m 4 cr 4 ) = (-)™ 1 " 



2C 



x S, 



J2 F 



k=0 



{k) (2£+ll 

(2k + iy 



x {£Q£Q\kQ) 2 (£-m l £m 3 \km 3 ~m l ) 

x (£ — m 2 £ m^k 777.4— m 2 ) . (A. 6) 

The fact that CG (£ £ 0\k 0) vanishes unless 21 + k is 
even and k preserves the well-known triangle condition, 
namely < k < 2£, enables to restrain summation to 
even values of k from to 2£, i.e. limited to k = 6 for / 
states, requiring only Slater integrals to F^ 6 ' , given 
in Table I for Plutonium. Replacement of result (| A.6|) 
into expression (|A.2[) , associated with selection rule for 
CG coefficients (£ m a\j M) which vanish unless M — 
777 + a, and final use of identity for ±| spins, (— )°'i~ <T2 = 
(2(7i) (2^2), one derives directly expression: 

^1234 — 0M 1 +M 2 ,M3+M 4 { — ) 

x 5^(2<7i)(2o- 2 )^Mi-tri<r 1 |jiMi) 



x M 2 -a 2 a 2 (is M 2 )(£ M 3 -a x a x \j 3 M 3 ) 

X ^M 4 -0-2(72b'4M4) 

X -(Mi -at) £, M3-CT1 1 fc M3-M1 } 

X (i-(M2-CT 2 )£M 4 -(72|fcM4-M2) , (A.7) 



from which one obtains expression (|2.1ip when we make 
an identification of definition (|2.12p for auxiliary function 
B i(jm;j'm'). 

Application of formula (|2.11[) to direct and exchange 
cases lead us to arrive at Eqs. (|2.16[) and (|2.17[) . re- 
spectively From these expressions and the CG related 
symmetries properties of function B^(jm; j'm') given 
in (|2.12[) . one can give the symmetry relations that a 
and &'s obey, reducing the number entries of Tables II 
and III. The matrix elements not given in those Tables 
can be consequently obtained using the following sym- 
metries: 

40M m l;j'M m 'l) = a{(jm;j'm') , (A.8) 



b{(j-m;j',m') = b{{j,m;j' -m') . (A.9) 



TABLE II: Values of a{{jm; f m') for the direct 
interaction matrix element {jmj'm'\V\jmj'm') = 
IZfc a k(j m ',j' m ')F f° r ^ — 3 in the jrrij basis. 



7 
J 


1 77X V 1 j' 

7 1 J 


\tyia 1 

1 7 




F W 


F (e) 


7/2 


7/2 7/2 


7/2 


49/441 49/5929 25/184041 




7/2 


5/2 


7 


-91 


-125 




7/2 


3/2 


-21 


-21 


225 




7/2 


1/2 

/. 


-35 


63 


-125 




5/2 


5/2 


1 


169 


625 




5/2 


3/2 


-3 


39 


-1125 




5/2 


1/2 


-5 


-117 


625 




3/2 


3/2 


9 


o 
y 


zUzo 




3/2 


1/2 


15 


-27 


-1125 




1/2 


1/2 


25 


81 


625 


7/2 


7/2 5/2 


5/2 


70/735 


7/1617 







7/2 


3/2 


-14 


-21 







7/2 


1/2 


-56 


14 







5/2 


5/2 


10 


-13 







5/2 


3/2 


-2 


39 







5/2 


1/2 


-8 


-20 







3/2 


5/2 


-30 


-3 







3/2 


3/2 


6 


9 







3/2 


1/2 


24 


-6 







1/2 


5/2 


-50 


9 







1/2 


3/2 


10 


-27 







1/2 


1/2 


40 


54 





5/2 


5/2 5/2 


5/2 100/ 1225 


1/441 







5/2 


3/2 


-20 


-3 







5/2 


1/2 


-80 


2 







3/2 


3/2 


4 


9 







3/2 


1/2 


16 


-6 







1/2 


1/2 


64 


4 






TABLE III: Values of b{{jm;j'm') for the ex- 
change interaction matrix element {jmj'm'\ V\j'm' jm) — 
Efc b{{jm-j'm')F {k) for I = 3 in the jm, basis. 



7/2 ±7/2 7/2 ±7/2 49/441 49/5929 25/184041 

±7/2 ±5/2 42 140 150 

±7/2 ±3/2 14 210 500 

±7/2 ±1/2 196 1200 



10 



±7/2 


±1/2 





98 


2250 


±7/2 


=F3/2 








3300 


±7/2 


±5/2 








3300 


±7/2 


±7/2 











±5/2 


±5/2 


1 


169 


625 


±5/2 


±3/2 


32 


60 


1400 


±5/2 


±1/2 


30 


2 


2100 


±5/2 


±1/2 





112 


2100 


±5/2 


±3/2 





210 


1050 


±5/2 


±5/2 











±3/2 


±3/2 


9 


9 


2025 


±3/2 


±1/2 


GO 


108 


1750 


±3/2 


±1/2 


40 


96 


700 


±3/2 


±3/2 











±1/2 


±1/2 


25 


81 


625 


±1/2 


±1/2 











7/2 ±7/2 5/2 ±5/2 175/11025 210/53361 25/184041 


±7/2 


±3/2 


140 


756 


200 


±7/2 


±1/2 





1323 


900 


±7/2 


±1/2 





1176 


3000 


±7/2 


±3/2 








8250 


±7/2 


±5/2 








19800 


±5/2 


±5/2 


150 


600 


150 


±5/2 


±3/2 


5 


1014 


875 


±5/2 


±1/2 


160 


486 


2800 


±5/2 


±1/2 





21 


6300 


±5/2 


±3/2 





1344 


10500 


±5/2 


±5/2 








11550 


±3/2 


±5/2 


75 


1000 


525 


±3/2 


±3/2 


90 


640 


2250 


±3/2 


±1/2 


30 


1 


5250 


±3/2 


±1/2 


120 


578 


8400 


±3/2 


±3/2 





686 


9450 


±5/2 


±5/2 





560 


6300 


±1/2 


±5/2 


20 


1200 


1400 


±1/2 


±3/2 


121 


120 


4375 


±1/2 


±1/2 


12 


300 


7500 


±1/2 


±1/2 


98 


735 


8750 


±1/2 


±3/2 


64 


60 


7000 


±1/2 


±5/2 





1050 


3150 


5/2 ±5/2 5/2 ±5/2 
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±3/2 
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4 
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±1/2 


60 
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±1/2 





14 
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±3/2 





14 





±5/2 


±5/2 











±3/2 


±3/2 


4 


9 





±3/2 


±1/2 


48 


10 





±3/2 


±1/2 


108 


5 





±3/2 


±3/2 











±1/2 


±1/2 


64 
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±1/2 


±1/2 












